#include "fuel_cell_infos.h"

#include <stdio.h>
#include <math.h>
#include "point.h"
//~ #include "mpi_tools.h"
#include "matrix.h"
#include "physical_constants.h"
#include "params.h"
#include "list.h"

#define PHYS_X1 0
#define PHYS_X2 0
#define PHYS_X3 0
#define PHYS_X4 0

tData MembraneQm(tData Tm){
    return PHYS_I * PHYS_I * sigma_m(Tm);
}

tData sigma_m(tData Tm){
    return (PHYS_X1 * lambda_m() - PHYS_X2) * exp( PHYS_X3 * (1/PHYS_X4 - 1/Tm));
}

tData lambda_m(float Ta, float Tc){
    tData a_a = (SH2O / MH2O - PHYS_R * Ta) / (A * exp( - B / Ta) );

    tData lambdaA;

    if( a_a >= 3 )
        lambdaA = 16.8;
    else
        if( a_a > 1 )
            lambdaA = 14 + 1.4 * (a_a -1);
        else
            lambdaA = 0.043 + 17.81 * a_a - 39.85 * a_a * a_a + 36 * a_a * a_a * a_a;

    tData a_c = (SH2O / MH2O) * PHYS_R * Tc / (A * exp(-B / Tc));

    if( a_c >= 3 )
        lambdaC = 16.8;
    else
        if( a_c > 1 )
            lambdaC = 14 + 1.4 * (a_c -1);
        else
            lambdaC = 0.043 + 17.81 * a_c - 39.85 * a_c * a_c + 36 * a_c * a_c * a_c;

    tData Dw;

    if( lambdaC < 2 )
        Dw = 0.000001;
    else
        if( lambdaC <= 3 )
            Dw = 0.000001 * ( 1 + 2 * (lamdaC - 3));
        else
            if( lamdaC <= 4.5 )
                Dw = 0.000001 ( 3 - 1.67 * (lambdaC - 3));
            else
                Dw = 1.25 * 0.000001;

        
}

